(* Mathematica package *)

Options[InputTransferC] = {
				LocalFieldsFromLimits -> False,
				SheetThickness -> 0, 
				InternalTransferCoefficientMethod -> "SII"(*"PartialSystemTransferCoefficients"*), 
				Incidence -> First, 
				InterfaceSide -> First, 
				AzimuthAngle -> 0};

InputTransferC[{Ep_,Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_, opts:OptionsPattern[]]:=
	
	Block[{ds, thetas, ns, inSI, inSII, inPhiv, inTin, 
		   voffset = Boole[!OptionValue[LocalFieldsFromLimits]] * OptionValue[InterfaceSide] @ {-1, 0},
		   (*TCM = OptionValue[InternalTransferCoefficientMethod],*)
		   SysSide = OptionValue[Incidence] @ {"0-", "kp1+'"},
		   vPos = OptionValue[InterfaceSide] @ {{"v+", "v-"}, {"v+'", "v-'"}}},
		{ds, thetas, ns} = 
			If[OptionValue[LocalFieldsFromLimits],
				{Insert[Dlist, OptionValue[SheetThickness], v], Insert[ThetaLIST, RefractionAngle[ThetaLIST[[1]], nlist[[1]], nbarINT], v+1], Insert[nlist, nbarINT, v+1]},
				{Dlist, ThetaLIST, nlist}];
		({inSI[#], inSII[#]}=
		 {AbelesSIv[#, Omega, thetas, ns, ds, v+voffset], AbelesSIIv[#, Omega, thetas, ns, ds, v+voffset]})& /@{"s", "p"};

		inPhiv = PhaseMv[Omega, thetas[[v+1+voffset]], ns[[v+1+voffset]], Quiet @ ds[[v+voffset]]];

		Table[inTin[vp] = PadRight[Evaluate[InternalTransferC[SysSide, vp, inSI[#], inPhiv, inSII[#],InternalTransferCoefficientMethod->OptionValue[InternalTransferCoefficientMethod]]& /@{"p", "s"}], 3], {vp, vPos}];
		

		(*Zrot[OptionValue[AzimuthAngle]] . *)J[ns[[v+1+voffset]], nbarINT] . (Pr[-thetas[[v+1+voffset]]] . inTin[vPos[[1]]] + Pr[thetas[[v+1+voffset]]] . inTin[vPos[[2]]]) * {Ep, Es, Ep}];

InputTransferC[{Ep_,Es_}, Omega_, Theta0_, n0_, ND:{{_, _}..}, nkp1_, v_Integer, nbarINT_, opts:OptionsPattern[]]:=
	Block[{ThetaLIST = {Theta0}~Join~RefractionAngle[Theta0,n0,ND[[All,1]]~Join~{nkp1}]},
		InputTransferC[{Ep,Es}, Omega, ThetaLIST, {n0}~Join~ND[[All,1]]~Join~{nkp1}, ND[[All,2]], v, nbarINT, opts]];

InputTransferC[{Ep_,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _}.. {_, _}}, v_Integer, nbarINT_, opts:OptionsPattern[]]:=
	InputTransferC[{Ep,Es}, Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]], v, nbarINT, opts];


InputTransferC[___] /; Message[General::badargs, InputTransferC] :=  "unevaluated";


InputTransferC[{0, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] :=
	Plus @@ (InternalTransferC["0-", "vTopSide", 
	   			   AbelesSIv["s", Omega, ThetaLIST, nlist, Dlist, v-1], 
		 	  	   PhaseMv[Omega, ThetaLIST[[v]], nlist[[v]], Dlist[[v-1]]], 
		 	   	   AbelesSIIv["s", Omega, ThetaLIST, nlist, Dlist, v-1]]) * {0,Es,0};

InputTransferC[{0, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] /; v == 1 :=
	(ReflectTopS[AbelesS["s", Omega, ThetaLIST, nlist, Dlist]] + 1) * {0, Es, 0};

InputTransferC[{0, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] /; v > Length @ Dlist :=
	TransmitTopS[AbelesS["s", Omega, ThetaLIST, nlist, Dlist]] * {0, Es, 0};


InputTransferC[{Ep_, 0}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] :=
	Ep * {Cos[ThetaLIST[[v]]] * Plus @@ #, 
		  0, 
		  (nlist[[v]] / nbarINT)^2 * Sin[ThetaLIST[[v]]] * Subtract @@ #}&[InternalTransferC["0-","vTopSide", 
																		   		 AbelesSIv["p",Omega,ThetaLIST,nlist,Dlist,v-1], 
																		   		 PhaseMv[Omega,ThetaLIST[[v]],nlist[[v]],Dlist[[v-1]]], 
																		   		 AbelesSIIv["p",Omega,ThetaLIST,nlist,Dlist,v-1]]];

InputTransferC[{Ep_, 0}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] /; v == 1 :=
	Ep * {Cos[ThetaLIST[[v]]] * (1 + #), 
		  0, 
		  (nlist[[v]] / nbarINT)^2 * Sin[ThetaLIST[[v]]] * (1-#)}& [ReflectTopS[AbelesS["p", Omega, ThetaLIST, nlist, Dlist]]];

InputTransferC[{Ep_, 0}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] /; v > Length @ Dlist :=
	Ep * TransmitTopS[AbelesS["p", Omega, ThetaLIST, nlist, Dlist]] * {Cos[ThetaLIST[[v+1]]], 0, (nlist[[v+1]] / nbarINT)^2 * Sin[ThetaLIST[[v+1]]]};
	
	
InputTransferC[{Ep_, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, v_, nbarINT_] :=
	InputTransferC[{Ep, 0}, Omega, ThetaLIST, nlist, Dlist, v, nbarINT] + InputTransferC[{0, Es}, Omega, ThetaLIST, nlist, Dlist, v, nbarINT];


(*Optimization for arguments of the form:
InputTransferC[{Ep_,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_]*)

InputTransferC[{0,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _}.., {_, _}}, v_Integer, nbarINT_]:=
	Plus @@ (InternalTransferC["0-", "vTopSide", 
	   			   AbelesSIv["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]], v-1], 
		 	  	   PhaseMv[Omega, ThetaND[[v,1]], ThetaND[[v,2]], ThetaND[[v, 3]]], 
		 	   	   AbelesSIIv["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]], v-1]]) * {0,Es,0};
	
InputTransferC[{0,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_] /; v == 1 :=
	(ReflectTopS[AbelesS["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]]]] + 1) * {0, Es, 0};

InputTransferC[{0,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_] /; v > Length @ ThetaND - 2 :=
	TransmitTopS[AbelesS["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]]]] * {0, Es, 0};
	


InputTransferC[{Ep_,0}, Omega_, ThetaND:{{_, _}, {_, _, _}.. {_, _}}, v_Integer, nbarINT_] :=
	Ep * {Cos[ThetaND[[v,1]]] * Plus @@ #, 
		  0, 
		  (ThetaND[[v,2]] / nbarINT)^2 * Sin[ThetaND[[v,1]]] * Subtract @@ #}&[InternalTransferC["0-","vTopSide", 
																		   		 AbelesSIv["p",Omega,ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]],v-1], 
																		   		 PhaseMv[Omega,ThetaND[[v,1]],ThetaND[[v,2]],ThetaND[[v,3]]], 
																		   		 AbelesSIIv["p",Omega,ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]],v-1]]];

InputTransferC[{Ep_,0}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_] /; v == 1 :=
	Ep * {Cos[ThetaND[[v, 1]]] * (1 + #), 
		  0, 
		  (ThetaND[[v,2]] / nbarINT)^2 * Sin[ThetaND[[v, 1]]] * (1-#)}& [ReflectTopS[AbelesS["p", Omega, ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]]]]];

InputTransferC[{Ep_,0}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_] /; v > Length @ ThetaND - 2 :=
	Ep * TransmitTopS[AbelesS["p", Omega, ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]]]] * {Cos[ThetaND[[v+1,1]]], 0, (ThetaND[[v+1,2]] / nbarINT)^2 * Sin[ThetaND[[v+1,1]]]};
	
	
InputTransferC[{Ep_, Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_] :=
	InputTransferC[{Ep,0}, Omega, ThetaND, v, nbarINT] + InputTransferC[{0,Es}, Omega, ThetaND, v, nbarINT];



(*******************************************************************************************************)
(*******************************************************************************************************)
(*******************************************************************************************************)
(*******************************************************************************************************)
(*Switch order of nbarINT and v args for convenience in GenerateInterfacialModel*)


(*InputTransferC[{Ep_,Es_}, Omega_, Theta0_, n0_, ND:{{_, _}..}, nkp1_, nbarINT_, v_Integer, opts:OptionsPattern[]]:=
	Block[{ThetaLIST = {Theta0}~Join~RefractionAngle[Theta0,n0,ND[[All,1]]~Join~{nkp1}]},
		InputTransferC[{Ep,Es}, Omega, ThetaLIST, {n0}~Join~ND[[All,1]]~Join~{nkp1}, ND[[All,2]], v, nbarINT, opts]];

InputTransferC[{Ep_,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer, opts:OptionsPattern[]]:=
	InputTransferC[{Ep,Es}, Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]], v, nbarINT, opts];

(*InputTransferC[{Ep_,Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer, opts:OptionsPattern[]]:=
	Block[{ds, thetas, ns, inSI, inSII, inPhiv, inTin, 
		   voffset = Boole[!OptionValue[LocalFieldsFromLimits]] * OptionValue[InterfaceSide] @ {-1, 0},
		   TCM = OptionValue[InternalTransferCoefficientMethod],
		   SysSide = OptionValue[Incidence] @ {"0-", "kp1+'"},
		   vPos = OptionValue[InterfaceSide] @ {{"v+", "v-"}, {"v+'", "v-'"}}},
		{ds, thetas, ns} = 
			If[OptionValue[LocalFieldsFromLimits],
				{Insert[Dlist, OptionValue[SheetThickness], v], Insert[ThetaLIST, RefractionAngle[ThetaLIST[[1]], nlist[[1]], nbarINT], v+1], Insert[nlist, nbarINT, v+1]},
				{Dlist, ThetaLIST, nlist}];
		({inSI[#], inSII[#]}=
		 {AbelesSIv[#, Omega, thetas, ns, ds, v+voffset], AbelesSIIv[#, Omega, thetas, ns, ds, v+voffset]})& /@{"s", "p"};

		inPhiv = PhaseMv[Omega, thetas[[v+1+voffset]], ns[[v+1+voffset]], Quiet @ ds[[v+voffset]]];

		Table[inTin[vp] = PadRight[Evaluate[tauIN[TCM, SysSide, vp, inSI[#], inPhiv, inSII[#]]& /@{"p", "s"}], 3], {vp, vPos}];

		(*Zrot[OptionValue[AzimuthAngle]] . *)J[ns[[v+1+voffset]], nbarINT] . (Pr[-thetas[[v+1+voffset]]] . inTin[vPos[[1]]] + Pr[thetas[[v+1+voffset]]] . inTin[vPos[[2]]]) * {Ep, Es, Ep}];*)


InputTransferC[{0, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] :=
	Plus @@ (InternalTransferC["0-", "vTopSide", 
	   			   AbelesSIv["s", Omega, ThetaLIST, nlist, Dlist, v-1], 
		 	  	   PhaseMv[Omega, ThetaLIST[[v]], nlist[[v]], Dlist[[v-1]]], 
		 	   	   AbelesSIIv["s", Omega, ThetaLIST, nlist, Dlist, v-1]]) * {0,Es,0};

InputTransferC[{0, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] /; v == 1 :=
	(ReflectTopS[AbelesS["s", Omega, ThetaLIST, nlist, Dlist]] + 1) * {0, Es, 0};

InputTransferC[{0, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] /; v > Length @ Dlist :=
	TransmitTopS[AbelesS["s", Omega, ThetaLIST, nlist, Dlist]] * {0, Es, 0};


InputTransferC[{Ep_, 0}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] :=
	Ep * {Cos[ThetaLIST[[v]]] * Plus @@ #, 
		  0, 
		  (nlist[[v]] / nbarINT)^2 * Sin[ThetaLIST[[v]]] * Subtract @@ #}&[InternalTransferC["0-","vTopSide", 
																		   		 AbelesSIv["p",Omega,ThetaLIST,nlist,Dlist,v-1], 
																		   		 PhaseMv[Omega,ThetaLIST[[v]],nlist[[v]],Dlist[[v-1]]], 
																		   		 AbelesSIIv["p",Omega,ThetaLIST,nlist,Dlist,v-1]]];

InputTransferC[{Ep_, 0}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] /; v == 1 :=
	Ep * {Cos[ThetaLIST[[v]]] * (1 + #), 
		  0, 
		  (nlist[[v]] / nbarINT)^2 * Sin[ThetaLIST[[v]]] * (1-#)}& [ReflectTopS[AbelesS["p", Omega, ThetaLIST, nlist, Dlist]]];

InputTransferC[{Ep_, 0}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] /; v > Length @ Dlist :=
	Ep * TransmitTopS[AbelesS["p", Omega, ThetaLIST, nlist, Dlist]] * {Cos[ThetaLIST[[v+1]]], 0, (nlist[[v+1]] / nbarINT)^2 * Sin[ThetaLIST[[v+1]]]};
	
	
InputTransferC[{Ep_, Es_}, Omega_, ThetaLIST_, nlist_, Dlist_, nbarINT_, v_Integer] :=
	InputTransferC[{Ep, 0}, Omega, ThetaLIST, nlist, Dlist, nbarINT, v] + InputTransferC[{0, Es}, Omega, ThetaLIST, nlist, Dlist, nbarINT, v];


(*Optimization for arguments of the form:
InputTransferC[{Ep_,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, v_Integer, nbarINT_]*)

InputTransferC[{0,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer]:=
	Plus @@ (InternalTransferC["0-", "vTopSide", 
	   			   AbelesSIv["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]], v-1], 
		 	  	   PhaseMv[Omega, ThetaND[[v,1]], ThetaND[[v,2]], ThetaND[[v, 3]]], 
		 	   	   AbelesSIIv["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]], v-1]]) * {0,Es,0};
	
InputTransferC[{0,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer] /; v == 1 :=
	(ReflectTopS[AbelesS["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]]]] + 1) * {0, Es, 0};

InputTransferC[{0,Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer] /; v > Length @ ThetaND - 2 :=
	TransmitTopS[AbelesS["s", Omega, ThetaND[[All,1]], ThetaND[[All,2]], ThetaND[[2;;-2,3]]]] * {0, Es, 0};
	


InputTransferC[{Ep_,0}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer] :=
	Ep * {Cos[ThetaND[[v,1]]] * Plus @@ #, 
		  0, 
		  (ThetaND[[v,2]] / nbarINT)^2 * Sin[ThetaND[[v,1]]] * Subtract @@ #}&[InternalTransferC["0-","vTopSide", 
																		   		 AbelesSIv["p",Omega,ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]],v-1], 
																		   		 PhaseMv[Omega,ThetaND[[v,1]],ThetaND[[v,2]],ThetaND[[v,3]]], 
																		   		 AbelesSIIv["p",Omega,ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]],v-1]]];

InputTransferC[{Ep_,0}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer] /; v == 1 :=
	Ep * {Cos[ThetaND[[v, 1]]] * (1 + #), 
		  0, 
		  (ThetaND[[v,2]] / nbarINT)^2 * Sin[ThetaND[[v, 1]]] * (1-#)}& [ReflectTopS[AbelesS["p", Omega, ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]]]]];

InputTransferC[{Ep_,0}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer] /; v > Length @ ThetaND - 2 :=
	Ep * TransmitTopS[AbelesS["p", Omega, ThetaND[[All,1]],ThetaND[[All,2]],ThetaND[[2;;-2,3]]]] * {Cos[ThetaND[[v+1,1]]], 0, (ThetaND[[v+1,2]] / nbarINT)^2 * Sin[ThetaND[[v+1,1]]]};
	
	
InputTransferC[{Ep_, Es_}, Omega_, ThetaND:{{_, _}, {_, _, _} .., {_, _}}, nbarINT_, v_Integer] :=
	InputTransferC[{Ep,0}, Omega, ThetaND, v, nbarINT] + InputTransferC[{0,Es}, Omega, ThetaND, v, nbarINT];
	
InputTransferC[___] /; Message[General::badargs, InputTransferC] :=  "unevaluated";

	
*)